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INTRODUCTION 

Cluster  orbits  are  defined  as  the  relative  trajectories  of  objects  traveling  through  space  in  close  proximity  to  each 
other.  Clusters  occur  when  satellites  in  similar  orbits  approach  each  other,  groups  of  satellites  fly  in  formation  to 
perform  specific  mission  functions,  or  several  objects  are  launched  or  deployed  from  the  same  object.  Each  case  has 
its  own  unique  nuances,  but  the  dynamics  of  how  the  objects  move  relative  to  each  other  are  functionally  the  same. 

The  number  of  satellites  in  the  geosynchronous  belt  has  been  steadily  increasing  for  almost  40  years.  While  the 
number  of  satellites  has  been  increasing,  the  space  in  the  geosynchronous  belt  has  not.  This  has  led  to  more  clusters 
of  satellites  operating  in  close  proximity  of  each  other  and  often  unintentionally  passing  within  kilometers  of  their 
neighbors'.  Clusters  of  objects  at  this  great  altitude  have  been  a  long-standing  challenge  to  space  surveillance  where 
the  satellites  can  be  mistaken  for  each  other. 

In  recent  years,  there  has  been  increasing  interest  in  the  use  of  satellites  flying  in  formation.  Several  missions  and 
mission  statements  have  identified  formation  flying  as  a  means  of  reducing  cost  and  adding  flexibility  to  space  based 
programs  or  to  accomplish  goals  that  are  not  possible  or  very  difficult  to  accomplish  with  a  single  satellite.  These 
missions  include  NASA's  Earth  Observing- 1  flying  in  formation  with  Landsat-7  and  several  European  missions.  In 
addition  NASA  has  developed  over  20  concepts  for  future  missions  involving  formation  flying.  Many  of  these 
missions  involve  highly  eccentric  orbits. 

During  ballistic  missile  launches,  several  objects  can  achieve  low  earth  orbit  including  the  reentry  vehicle,  booster, 
fairings,  and  the  possibility  of  several  balloon-like  decoys  attempting  to  confuse  missile  defense  systems.  It  is  the 
responsibility  of  the  missile  defense  system  to  track  the  cluster  of  objects  and  quickly  identify  the  reentry  vehicle. 

The  Cluster  Orbits  with  Perturbations  of  Keplerian  Elements  (COWPOKE)  equations  can  provide  the  theoretical 
foundation  for  analysis  tools  supporting  a  variety  of  applications.  Optical  space  surveillance  often  reveals  multiple 
satellites  in  a  single  telescope  field  of  view.  An  understanding  of  the  relative  dynamics  may  allow  for  satellite 
identification  based  solely  on  relative  position.  Eurther  research  may  apply  these  results  to  relative  orbit 
determination  of  objects  in  a  cluster  or  during  close  approach  encounters.  This  approach  may  prove  to  be  more 
accurate  than  performing  orbit  determination  on  the  individual  objects  and  determining  differences.  The  equations 
might  also  be  used  to  support  the  missile  tracking  experiments.  If  tracking  sensors  have  limited  fields  of  view,  they 
may  be  required  to  track  ballistic  clusters  individually;  meaningful  relative  equations  of  motion  may  simplify  the 
transition  from  object  to  object.  Even  if  objects  are  tracked  simultaneously,  the  relative  equations  of  motion  can  be 
used  as  part  of  the  discrimination  process  to  identify  balloon-like  decoys.  The  COWPOKE  equations  would  also  be 
extremely  relevant  to  satellite  formation  flying  work,  such  as  the  initial  Air  Force  Research  Laboratory’s  TechSat  21 
program,  for  formation  flying  design,  analysis,  and  guidance  and  control  applications.  While  a  great  deal  of  work 
has  already  been  completed  in  support  of  formation  flying  missions,  none  provide  the  level  of  intuition  and 
understanding  as  one  that  uses  Keplerian  elements.  In  addition,  none  of  the  work  is  accurate  for  long-term  relative 
motion  for  orbits  with  eccentricities  up  to  0.7. 

A  simple  set  of  equations  describing  the  relative  motion  of  spacecraft  clusters  is  needed  for  analysis,  design,  and/or 
tracking  of  clusters.  Much  of  the  previous  work  in  this  field  has  relied  on  using  Hill's  equations^  (also  known  as  the 
Clohessy-Wiltshire  equations^).  Hill's  equations  describe  the  relative  motion  of  spacecraft  using  a  spacecraft- 
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centered  coordinate  system.  However,  Hill's  equations  assume  that  the  reference  orbit  is  circular,  the  objects  are 
close  together,  and  there  are  no  perturbations  to  simple  two-body  motion.  Several  researchers  have  pointed  out  the 
severe  limitations  inherit  in  these  assumptions.  To  rectify  this,  Gim  and  Alfriend*  used  energy  methods  to  develop  a 
state  transition  matrix  to  describe  relative  motion  in  non-circular  orbits  under  the  influence  of  perturbations.  While 
closed  form  analytic  solutions  in  a  transformed  variable  space  are  valuable  for  analysis,  these  approaches  may  be 
awkward  for  many  applications  since  transformations  are  required  to  go  from  the  canonical  variable  space  to  more 
traditional  representations  of  satellite  orbits  such  as  Keplerian  elements.  Additionally,  most  working  level  engineers 
lack  the  background  to  implement  this  approach.  Previous  work  by  Garrison  et  al.^  developed  equations  of  motion 
for  elliptical  orbits  in  terms  of  the  true  anomaly  instead  of  time.  Also,  these  equations  were  developed  with 
rendezvous  in  mind  and  are  developed  only  to  second  order  in  eccentricity.  In  addition,  Melton*’  developed  a  state 
transition  matrix  for  relative  motion  in  eccentric  orbits  that  is  time  dependent.  Melton  relied  on  some  of  the  same 
foundations  that  are  used  in  this  paper,  but  the  development  was  only  to  second  order  in  eccentricity  and  is  not 
accurate  for  orbits  of  high  eccentricity.  Baoyin  et  al,'^  did  similar  work  but  assumed  matching  orbital  periods;  this 
work  also  provides  physical  insight  into  formation  flying  in  near  circular  orbits.  This  paper  develops  physically 
meaningful  equations  of  relative  motion  for  space  objects  in  non-circular  orbits  using  Keplerian  elements. 

Previous  work  by  the  authors^  led  to  the  development  of  a  set  of  equations  that  describe  the  first  order  effects  of 
Earth  oblateness  on  the  relative  motion  of  objects  in  circular  orbits.  In  addition,  they  developed  a  simple  set  of 
equations  to  describe  the  effects  of  Earth  oblateness  for  polar  orbits.  Then,  they  further  developed  the  equations  to 
describe  the  motion  for  all  inclinations*.  Finally,  they  examined  the  long-term  evolution  of  the  relative  motion  for 
circular  orbits  and  presented  an  approach  to  describe  relative  motion  of  satellites  in  elliptical  orbits  without 
perturbations  and  assuming  matching  periods^.  This  last  step  provided  the  building  blocks  for  the  COWPOKE 
equations. 

In  the  past,  equations  of  relative  motion  such  as  Hill’s  equations  were  developed  by  differencing  the  quasi-inertial 
equations  of  motion  and  mapping  those  differences  into  the  rotating  coordinate  system.  The  result  was  a  set  of 
equations  that  describe  the  relative  differences  between  the  satellites  typically  in  terms  of  radial,  cross-track,  and 
along-track  components.  Rather  than  using  the  traditional  algebraic  approach,  the  COWPOKE  equations  are 
developed  using  a  geometric  approach.  Here,  the  geometric  properties  and  definitions  of  Keplerian  orbital  elements 
will  be  used  directly  to  map  orbital  element  differences  into  the  radial,  cross-track,  and  along-track  relative  motion. 
In  addition  to  the  increased  intuitiveness  of  a  geometric  approach,  by  using  Keplerian  orbital  elements  and 
differences  in  those  orbital  elements,  we  can  take  advantage  of  existing  perturbation  model  development  and 
incorporate  meaningful  dynamics  into  the  relative  equations  of  motion  much  easier  than  with  the  algebraic 
approaches. 

There  are  challenges  to  this  approach.  First  and  foremost,  the  geometric  properties  of  Keplerian  elements  are  a 
function  of  the  true  anomaly,  which  has  non-uniform  variation  with  time.  Ideally,  analytical  solutions  are  expressed 
in  terms  of  the  mean  anomaly,  which  does  vary  linearly  with  time  (outside  the  influence  of  perturbations).  The 
relationship  between  the  mean  anomaly  and  eccentric  anomaly  is  given  by  the  well-  known  transcendental  Kepler’s 
equation.  This  research  uses  a  series  expansion  to  approximate  the  true  anomaly  as  a  function  of  time.  A  relation  of 
this  sort  was  successfully  employed  as  part  of  the  building  blocks  of  the  COWPOKE  equations. 

Using  that  one  relationship,  one  can  construct  the  relative  motion  of  formation  flying  satellites  by  reconstructing  the 
orbits  individually  based  on  their  Keplerian  elements  and  then  differencing  the  two.  This  paper  develops  true 
equations  of  relative  motion  for  clusters  of  objects.  And,  unlike  the  previous  work,  this  research  incorporates 
perturbation  models  and  be  generalized  to  account  for  small  differences  in  semimajor  axis.  This  results  in  a  set  of 
equations  that  describe  the  radial,  cross-track,  and  along-track  differences  between  two  space  objects  in  formation  or 
a  cluster  based  on  Keplerian  orbital  element  differences. 
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METHODS,  ASSUMPTIONS,  AND  PROCEDURES 

The  formulation  of  the  COWPOKE  equations  has  three  major  challenges.  The  first  is  choosing  a  suitable  reference 
frame  for  the  relative  equations  of  motion  and  expressing  that  frame  in  terms  of  Keplerian  elements  and  element 
differences.  The  second  is  representing  the  geometrical  true  anomaly  as  a  function  of  time.  The  final  challenge  is 
incorporating  relevant  perturbation  effects.  Only  the  first  two  challenges  are  addressed  in  this  work. 

Classical  relative  motion  approaches  use  a  Cartesian  reference  frame  to  describe  relative  motion.  The  components 
of  this  frame  are  radial,  cross-track,  and  along-track  differences.  A  limitation  of  this  approach  is  that  significant 
cross-mapping  between  components  occurs  when  the  separation  distance  between  the  satellites  is  not  small.  For 
instance,  if  two  satellites  travel  together  in  a  circular  orbit  but  are  separated  by  0.01  rad  in  true  anomaly,  the  relative 
position  of  one  satellite  with  respect  to  the  other  will  have  an  along-track  separation  as  anticipated  but  will  also  have 
a  radial  component  despite  the  fact  that  the  satellites  are  at  the  same  altitude  and  within  the  same  orbit.  Relative 
position  difference  components  that  are  artifacts  of  the  chosen  coordinate  frame,  such  as  the  one  outlined  in  the 
preceding  example,  can  greatly  complicate  the  analysis  of  the  relative  motion. 

A  second  coordinate  frame  choice  is  to  use  spherical  separations.  Here,  the  relative  position  difference  is  described 
by  an  altitude  difference  from  a  sphere  having  the  radius  of  the  reference  satellite  and  angular  components 
perpendicular  to  and  along  the  reference  satellite’s  direction  of  motion  projected  onto  the  sphere.  These 
components,  5r,  and  Sat,  are  illustrated  in  Figure  1.  The  spherical  reference  frame  is  well-suited  to  describing 
the  position  of  satellites  in  circular  or  near-circular  orbits;  however,  when  the  orbital  eccentricity  becomes  large,  the 
spherical  reference  frame  suffers  from  the  same  limitations  as  the  Cartesian  system. 

A  third  option  to  describe  the  relative  motion  of  satellites  in  elliptic  orbits  is  to  use  an  ellipsoidal  reference  system. 
This  is  similar  to  the  spherical  system  except  that  the  components  are  mapped  along  an  ellipsoid  rather  than  a 
sphere.  This  method  has  geometric  challenges  describing  altitude  variations  and  can  still  result  in  difficult  to 
understand  position  differences. 

For  this  formulation,  we  chose  the  spherical  reference  frame  since  it’s  slightly  better  than  Cartesian  for  describing 
the  orbital  motion  and  provides  the  physical  insight  we  desire.  Additionally,  it  is  a  simple  geometric  mapping 
between  spherical  and  Cartesian  given  the  quantities  described  in  this  work  if  Cartesian  coordinates  are  desired.  In 
fact,  the  desired  reference  frame  will  likely  be  a  function  of  application,  but  it  is  hoped  that  enough  information  is 
provided  here  so  these  equations  can  be  easily  reformulated  in  other  reference  frames. 
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The  spherical  coordinate  system  used  here  can  be  realized  using  position  and  velocity  vectors  of  the  two  satellites: 


008(7)=  ~ 


cos(y) 


P. 


8r  =  r^-r, 


□  □ 

V,  =v,- 


COS((|)): 


V|  T, 


V  '1  J' 


8/  Vj 

sin(p)  =  sin(7)sin((|)) 
cos(y) 


cos(x)  = 


cos(P) 


§r  =  -  fj 

duct  =  rjP 
8at=  rfii 


(1) 


Here,  the  auxiliary  vector  quantities,  denoted  by  primes,  are  introduced  so  right  planar  triangles  can  be  used  in  the 
determination  of  desired  angular  values.  The  spherical  components  can  be  expressed  in  terms  of  Keplerian  elements 
by  using  the  following  definitions'”: 


□ 

r  = 


2  cos(Q)cos^  +  /)  -sin(Q)sin(co+  /)cos(i) 
sin(fi)cos^+  /)-cos(P)sin(co+ /)cos(i) 

l  +  ecos(f) 

L  sin(C0  + /)sin(0 


(2) 


□ 

V  = 


-p 


cos(P)(sin(co  +  /)  +e  sin(co))+  sin(Q)(cos(co+  /)+  ecos^))cos(0 
sin(Q)(sin(C0  +  /)  +  e  sin(co))+  cos(P)(cos(co+  /)+  ecos^))cos(0 
— (cos(co  + /)+ ecos(co))  sin(0 


However,  the  resulting  equations  become  difficult  to  manipulate.  Instead,  a  simple  geometrical  mapping  of 
Keplerian  element  and  element  differences  into  the  spherical  components  of  relative  motion  can  be  achieved  with 
the  following: 


_  (a,  +5fl)(l-(e,  +5e)^)  aiCl-e,^) 

1+ (e, +5e)cos(/j  + 5/)  1+e,  cos(/',) 

=  P  =-6r2sin(t[)cos(C0j +5co+/i +50  +  5/8)11(0)1  +  50)+ /i  +  pC)  (3) 

=a  =  (5o)  + 500085/)+ 5Qco8(/i) 

/  ^ 
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The  spherical  cross-track  and  along-track  terms  are,  in  fact,  approximations  accurate  to  first  order  in  the  Keplerian 
elements  differences.  Figure  2  illustrates  how  these  terms  are  derived. 


Figure  2:  Spherical  Components  in  Terms  of  Keplerian  Elements  and  Element  Differences 


We  can  further  simplify  the  radial  expression  to  first  order  in  the  Keplerian  element  differences: 

s  s  ^-fl(2ei+(l  +  ei^)cos(/  ))  a  e  (1-e  ^)sin(/ ) 

or= - i - oa-l - i i 5 oe-l-— ^ - T^of  (4) 

l  +  e,cos(/-,)  {l  +  e,cos(f,)f  (l  +  e,cos(f,))^ 

From  these  Keplerian  element  representations  of  the  relative  motion,  we  can  derive  some  physical  insight  into  how 
orbital  differences  result  in  relative  motion.  If  we  assume  a  circular  reference  orbit,  Eqs.  (3)  and  (4)  resemble  the 
solutions  to  Hill’s  equations  where  the  cross-track  motion  becomes  a  simple  oscillation,  the  altitude  difference  maps 
directly  into  the  radial  motion  with  a  periodic  variation  introduced  by  the  eccentricity,  and  the  altitude  difference 
couples  into  the  along-track  component  through  the  true  anomaly  resulting  in  a  combination  of  secular  and  periodic 
effects’. 

While  the  equations  above  are  useful  and  do  help  provide  physical  insight  into  the  relative  motion,  they  are 
expressed  in  terms  of  true  anomaly.  To  express  the  relative  motion  as  a  function  of  time,  one  must  express  the 
above  equations  in  terms  of  mean  anomaly.  At  a  fundamental  level,  this  entails  finding  a  suitable  approximate 
solution  to  Kepler’s  equation.  Battin  provides  several  methods  for  doing  so’”  For  the  initial  COWPOKE 
development,  it  was  decided  to  use  a  Fourier-Bessel  expansion  of  the  true  anomaly  in  terms  of  the  mean  anomaly 
and  eccentricity. 

Using  the  geometric  properties  of  the  true  and  eccentric  anomalies  and  Kepler’s  equation,  the  following  relationship 
can  be  found: 


df _ i/l-  e’ 

dM  (1-ecosi?)^ 
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This  can  be  expressed  in  terms  of  a  Fourier  cosine  series  that  introduces  the  Bessel  functions,  J„: 
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dM 


Bo+^B^cos(kM) 

V  fci  y 


Bo  = 
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B,= 
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{-ke)  =(-ke)" 


y  i-iri-kef" 

h2^’^''fn\{n  +  m)\  ’ 


^n(lfc  +  „l+2/-y) 


./=i 


(6) 


Integrating  this  equation  with  respect  to  the  mean  anomaly  provides  an  expression  for  the  true  anomaly  in  terms  of 
the  mean  anomaly: 


/=M  +  2£| 

A'=l  ^ 


£j„(-^e)|3 


|r+„| 


sm{kM) 


(V) 


The  result  of  this  expansion  is  a  sine  series  with  coefficients  that  are  power  series  in  eccentricity.  The  lowest  order 
of  eccentricity  for  each  coefficient  series  is  k;  thus  the  upper  limit  for  k  can  be  chosen  based  on  the  eccentricity  of 
the  reference  orbit  and  the  desired  accuracy.  Figure  3  plots  eccentricity  to  the  (k+1)  power.  The  figure  is  meant  to 
provide  an  indicator  for  the  relative  accuracy  of  Eq.  (7)  for  a  given  eccentricity  and  choice  of  k.  For  instance,  if  one 
had  a  references  eccentricity  of  0.3  and  desired  series  truncation  errors  below  1%,  then  one  would  choose  to  truncate 
the  series  at  k=4  or  the  fourth  power  of  eccentricity.  When  k=8,  one  sees  errors  below  10%  up  to  e=0.7. 
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Figure  3:  Relative  Error  Order  of  Magnitude  for  Eccentricity  Series 
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Evaluating  Eq.  (7)  up  to  the  ninth  frequency  term  (8M)  yields: 
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In  examining  Eq.  (8),  one  can  see  that  the  leading  term  in  each  of  the  eccentricity  series  does  not  decrease  for  all 
eccentricities;  however,  the  second  term  in  the  eccentricity  series  is  always  negative  and  serves  to  reduce  the 
magnitude  of  that  frequency’s  contribution  to  the  expansion  as  a  whole.  Thus,  the  amplitude  of  each  frequency  term 
decreases  as  the  frequency  increases  allowing  convergence  of  the  entire  expansion  if  enough  terms  are  included  in 
each  individual  eccentricity  series.  This  means  that  for  higher  values  of  eccentricity,  it  is  very  important  to  include 
at  least  two  or  more  terms  for  each  eccentricity  series.  In  Eq.  (8),  if  the  e’  and  e'°  terms  were  left  out  of  the  7M  and 
%M  frequency  terms,  the  error  due  to  this  truncation  could  be  larger  than  if  the  IM  and  HM  terms  were  not  included 
at  all.  For  higher  values  of  eccentricity,  care  must  be  taken  to  make  certain  that  enough  terms  are  included  in  each 
eccentricity  series  to  ensure  one  is  not  adding  error  with  a  given  frequency  term.  The  good  news  is  that  each 
eccentricity  series  converges  fairly  quickly  so  that  powers  of  eccentricity  at  lower  frequencies  are  much  less 
important  than  at  higher  frequencies;  this  means  that  one  can  get  away  with  fewer  terms  in  each  eccentricity  series 
even  if  higher  frequency  terms  are  required.  Additionally,  for  smaller  values  of  eccentricity,  0.3  and  below, 
truncation  becomes  much  less  of  an  issue  since  the  e^  factor  greatly  reduces  the  impact  of  the  higher  order  terms  in 
each  eccentricity  series  and  each  higher  frequency. 


Truncation  issues  aside,  the  Fourier-Bessel  series  expansion  of  the  true  anomaly  is  an  important  component  to  the 
formulation  of  the  COWPOKE  equations;  however,  an  additional  step  must  be  taken  to  map  mean  anomaly 
differences  into  true  anomaly  differences.  This  is  accomplished  by  simply  applying  a  first  order  perturbation  to  Eq. 
(8): 


(9) 

dM  de 

The  general  formulation  includes  higher  order  element  difference  terms  and  coupling  terms  between  the  element 
differences  but  for  cluster  orbits,  including  the  first  order  terms  should  provide  accuracy  to  two  orders  of  magnitude 
smaller  than  the  separation  distance.  Including  only  the  first  order  difference  terms  for  the  true  anomaly  expression 
yields; 
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The  truncation  issues  for  Eq.  (10)  are  similar  to  those  of  Eq.  (8)  but  are  further  compounded  by  the  exponential 
factors  produced  when  taking  the  partial  derivative  with  respect  to  eccentricity.  Note  that  the  e^  coefficient  is  twice 
large  as  the  e’  coefficient  in  the  8  A/  frequency  series  for  the  6e  component. 

For  unperturbed  satellite  motion,  we  can  express  the  mean  anomaly  as  a  function  of  time  5M  as  a  function  of  5M  at 
epoch,  semimajor  axis  difference,  and  time: 


Mi  =  Mo  + 
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If  one  chooses,  the  mean  motion  terms  in  the  above  expression  can  be  replaced  by  a  second  order  perturbation  of  the 
mean  motion  with  respect  to  &  without  significant  error  for  most  applications: 


5M  =  5Mg  +  Jj!  ai^^5a+  —  q 

V  2  8  J 
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Eqs.  (8),  (10),  and  (11)  can  now  be  substituted  into  Eq.  (3)  to  complete  the  COWPOKE  equations  for  unperturbed 
satellite  motion.  Here  are  the  COWPOKE  equations  with  only  first  order  eccentricity  terms  included  in  the 
expansion  of  the  true  anomaly  and  true  anomaly  difference: 
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Note  that  we  have  taken  the  linearized  form  of  the  radial  component  as  presented  in  Eq.  (4).  From  Figure  3,  we 
would  expect  that  these  equations  are  accurate  to  the  1  %  level  up  to  eccentricities  of  0. 1 . 


While  this  work  does  not  currently  include  perturbations  to  the  satellite  dynamics,  one  of  the  key  focuses  of  the 
COWPOKE  development  was  to  allow  for  the  inclusion  of  additional  force  model  effects.  This  can  be 
accomplished  by  expressing  the  Keplerian  elements  and  element  differences  as  functions  of  time.  This  will  be  the 
focus  of  future  work. 
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RESULTS  AND  DISCUSSION 


Simulations  were  performed  to  quantify  the  error  inherent  in  the  approximations  used  in  the  COWPOKE 
formulation  for  two  test  cases.  The  first  test  case  employed  a  near-circular,  low-Earth  orbit  (e=0.01)  and  the  second 
a  high  altitude  eccentric  orbit  (e=0.7).  The  test  cases  were  performed  in  the  MATLAB  environment  using  two-body 
dynamics;  it  is  understood  that  that  these  simulations  only  measure  the  effectiveness  of  the  COWPOKE  equations  to 
model  two-body  motion  and  are  not  indicative  of  the  performance  for  real  satellite  motion. 

Truth  data  for  the  simulations  were  determined  by  calculating  the  mean  anomaly  of  each  satellite  as  a  function  of 
time,  converting  theKeplerian  elements  to  position  and  velocity  vectors  (with  an  algorithm  taken  from  Vallado"), 
and  then  mapping  the  position  differences  into  the  spherical  radial,  cross-track,  and  along-track  components 
described  in  Eq.  (1).  The  COWPOKE  results  come  from  a  direct  mapping  of  the  Keplerian  element  and  element 
differences  into  those  components  using  the  methods  described  in  the  previous  section.  It  should  be  noted  that  the 
angular  cross-track  and  along-track  values,  a  and  p ,  were  used  for  comparison  purposes  in  the  simulations  rather 
than  the  arc-lengths,  &fand  &t.  Some  additional  error  in  the  arc-lengths  will  be  present  due  to  the  error  in 
estimating  the  magnitude  of  the  position  vector  of  the  reference  satellite,  but  that  is  typically  the  same  order  of 
magnitude  or  smaller  than  the  error  in  the  angular  separation. 

Table  1  contains  the  orbital  elements  and  element  differences  for  the  near-circular,  low-Earth  orbit  (LEO)  case;  the 
reference  satellite,  satellite  1,  has  the  initial  Keplerian  elements  given  in  the  “Reference  Elements”  column  while  the 
second  satellite,  satellite  2,  has  the  initial  Keplerian  elements  of  the  “Reference  Elements”  plus  the  “Element 
Differences.”  The  simulation  span  covers  2  hours  or  just  over  one  orbital  period.  The  COWPOKE  formulation 
includes  only  first  order  eccentricity  terms  and  uses  the  radial  component  approximation;  these  equations  were 
identical  toEq.  (13). 


Table  1 :  Keplerian  Elements  and  Element  Differences  for  LEO  Test  Case 


Reference 

Element 

Elements 

Differences 

a 

7000  km 

0.01  km 

e 

0.01 

0.01 

i 

0.785  rad  (45  deg) 

O.Olrad 

a. 

0  rad  (0  deg) 

0.01  rad 

03 

4.712rad(270deg) 

0.01  rad 

_ _ 

1.751  rad  (90  deg) 

0.01  rad 

Figures  4-6  plot  the  spherical  radial  difference,  angular  cross-track,  and  angular  along-track  separations, 
respectively.  Each  figure  actually  contains  two  plots:  the  first  shows  the  satellite  separations  described  by  the  truth 
orbits  and  the  COWPOKE  approximation,  and  the  second  shows  the  differences  between  the  truth  and  COWPOKE 
methods.  Thus,  the  second  plot  in  each  figure  is  the  error  inherent  in  the  given  COWPOKE  formulation. 

Figure  4  shows  the  COWPOKE  radial  error  to  be  around  the  2%  level.  Figures  5  and  6  show  the  COWPOKE  cross¬ 
track  and  along-track  errors  are  below  1%.  It  is  interesting  to  note  that  the  radial  and  cross-track  component  errors 
have  frequencies  that  appear  to  be  once  per  orbit  while  the  along-track  component  has  an  error  frequency  of  twice 
per  orbit.  We  test  the  COWPOKE  approximations  directly  to  determine  the  dominant  sources  of  these  errors. 

Before  continuing  with  the  error  analysis,  it  is  useful  to  examine  Figures  4-6  and  observe  the  physical  insight  into 
the  relative  motion  we  can  gain  from  the  COWPOKE  equations.  Eq.  (4)  tells  us  that  the  relative  motion  in  the  radial 
direction  is  dominated  by  a  once  per  orbit  signature  due  to  the  eccentricity  differences;  the  contributions  from  the 
semimajor  axis  difference  are  small,  and  the  contribution  from  the  true  anomaly  difference  is  an  order  of  eccentricity 
smaller  than  the  contribution  from  the  eccentricity  difference.  Eq.  (3)  shows  that  the  cross-track  motion  is 
comprised  of  once  per  orbit  signatures  from  the  inclination  and  right  ascension  of  the  ascending  node  differences. 
Eq.  (13)  indicates  that  the  along-track  motion  is  offset  from  zero  due  to  the  right  ascension,  argument  of  perigee,  and 
mean  anomaly  differences  but  also  has  a  once  per  orbit  periodic  signature  with  an  amplitude  of  twice  the  reference 
eccentricity. 
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Figure  4:  LEO  Test  Case  Radial  Differences  and  COWPOKE  Errors 
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Figure  5:  LEO  Test  Case  Angular  Cross-Track  Differences  and  COWPOKE  Errors 
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Figure  6:  LEO  Test  Case  Angular  Along-Track  Differences  and  COWPOKE  Errors 


Figure  7  plots  the  error  in  the  first  order  Fourier-Bessel  series  expansion  of  the  true  anomaly  of  the  reference 
satellite.  From  Figure  3,  we  expect  to  see  around  1%  error,  and  we  do.  Also  note  the  twice  per  orbit  signature  in  the 
error  plot;  this  is  also  expected  since  the  expansion  only  include  once  per  orbit  terms.  Figure  8  plots  the  error  in  the 
approximation  of  the  true  anomaly  difference;  the  error  in  this  approximation  is  similar  to  the  error  in  the  true 
anomaly  approximation.  One  sees  that  the  dominant  part  of  the  along-track  error  is  due  to  the  COWPOKE 
approximation  error  of  the  true  anomaly  difference.  If  one  chose,  this  error  could  be  reduced  by  adding  additional 
terms  to  the  series  approximation  for  this  variable. 


Figure  7:  True  Anomaly  Approximation  Error  for  LEO  Test  Case 
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Figure  8:  True  Anomaly  Difference  Approximation  Error  for  LEO  Test  Case 


Figure  9  plots  the  errors  in  the  COWPOKE  coordinate  system  approximations  for  the  LEO  test  case.  One  can  see 
that  the  dominant  part  of  the  radial  and  cross-track  errors  is  due  to  the  coordinate  system  approximations  made  in  the 
COWPOKE  formulations.  The  radial  component  error  can  be  reduced  by  using  the  analytical  formulation  for  the 
radial  difference  while  an  improved  description  of  the  cross-track  component  would  have  to  be  derived  to  reduce  the 
cross-track  error. 


Table  2  contains  the  orbital  elements  and  element  differences  for  the  highly  eccentric  Earth  orbit  (HEO)  test  case. 
The  simulation  spans  12  hours  or  just  over  one  orbital  period.  The  COWPOKE  formulation  includes  eighth  order 
eccentricity  up  to  frequency  6M  and  tenth  order  eccentricity  in  frequency  IM  and  8M;  these  are  all  of  the  terms 
included  in  Eqs.  (8)  and  (10).  This  test  case  does  not  use  the  radial  component  approximation. 

Table  2:  Keplerian  Elements  and  Element  Differences  for  HEO  Test  Case 


Reference 

Elements 

Element 

Differences 

a 

27000  km 

0.01  km 

e 

0.7 

0.01 

i 

0.785  rad  (45  deg) 

O.Olrad 

a. 

0  rad  (0  deg) 

0.01  rad 

0) 

4.712  rad  (270  deg) 

0.01  rad 

M„ 

1.751  rad  (90  deg) 

0.01  rad 

Figures  10-12  plot  the  spherical  radial  difference,  angular  cross-track,  and  angular  along-track  separations, 
respectively.  As  with  Figures  4-6,  the  first  plot  in  each  figure  shows  the  satellite  separations  described  by  the  truth 
orbits  and  the  COWPOKE  approximation,  and  the  second  plot  shows  the  differences  between  the  truth  and 
COWPOKE  methods.  Thus,  the  second  plot  in  each  figure  is  the  error  inherent  in  the  given  COWPOKE 
formulation  for  this  highly  eccentric  test  case. 
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Figure  9:  COWPOKE  Coordinate  System  Approximation  Errors  for  LEO  Test  Case 
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Figure  10:  HEO  Test  Case  Radial  Differences  and  COWPOKE  Errors 
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Figure  1 1 :  FIFO  Test  Case  Angular  Cross-Track  Differences  and  COWPOKE  Errors 
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Figure  12:  HEO  Test  Case  Angular  Along-Track  Differences  and  COWPOKE  Errors 


Eigure  10  shows  the  COWPOKE  radial  error  to  he  significant  with  maximum  errors  being  on  the  same  order  of 
magnitude  as  the  relative  altitude  differences.  Figures  1 1  and  12  show  the  COWPOKE  cross-track  and  along-track 
errors  are  on  the  order  of  10%.  The  sources  of  these  errors  can  be  determined  by  analyzing  the  approximations  used 
in  the  COWPOKE  formulation. 


Eigure  13  plots  the  error  in  the  Eourier-Bessel  series  expansion  of  the  true  anomaly  of  the  reference  satellite.  Erom 
Figure  3,  we  expect  to  see  around  10%  error,  and  we  do  along  with  a  mixture  of  frequencies  in  the  signature  due  to 
the  series  truncation.  The  error  in  the  true  anomaly  approximation  maps  into  the  frequency  terms  and  becomes  the 
dominant  source  of  error  for  the  cross-track  component.  Figure  14  plots  the  error  in  the  approximation  of  the  true 
anomaly  difference;  the  error  in  this  approximation  is  an  order  of  magnitude  smaller  than  the  error  in  the  true 
anomaly  approximation  since  we  have  removed  an  order  of  eccentricity  in  the  differentiation  process  and  replaced  it 
by  an  eccentricity  difference.  Since  the  eccentricity  difference  is  an  order  of  magnitude  smaller  than  the 
eccentricity,  the  eccentricity  difference  error  is  much  smaller.  One  can  see  that  the  dominant  part  of  the  radial  and 
along-track  errors  is  due  to  the  COWPOKE  approximation  error  of  the  true  anomaly  difference  and  likely  the 
eccentricity  series  truncation  for  the  8M  frequency  terms.  If  one  chose,  these  errors  could  be  reduced  by  adding 
higher  order  terms  to  the  series  approximations  for  these  variables. 
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Figure  13:  True  Anomaly  Approximation  Error  for  HEO  Test  Case 


Figure  14:  True  Anomaly  Difference  Approximation  Error  for  HEO  Test  Case 

Figure  15  plots  the  errors  in  the  COWPOKE  coordinate  system  approximations  for  the  HEO  test  case.  There  is  no 
plot  for  the  radial  component  since  the  radial  approximation  was  not  used,  and  thus  the  coordinate  system  error  is 
zero.  The  cross-track  and  along-track  components  remain  around  or  below  the  1%  level  even  for  the  high 
eccentricity  case.  It  would  seem  that  the  coordinate  system  gives  the  COWPOKE  equations  an  error  bound  of  about 
1%. 

The  cross-track  and  along-track  errors  for  the  HEO  case  are  in  line  with  expectations  given  the  series  truncations 
made  in  the  COWPOKE  formulation  for  this  test  case.  The  source  of  the  radial  error  was  larger  than  expected  or 
desired.  Eurther  analysis  of  Eq.  (4)  shows  that  there  is  coupling  between  the  true  anomaly  error  and  the  eccentricity 
difference.  The  only  way  to  mitigate  this  en'or  is  to  reduce  the  error  in  the  true  anomaly  difference  approximation 
which  can  be  accomplished  by  adding  additional  terms  to  the  series  expansion  or  by  solving  Kepler’s  equation 
iteratively.  This  would  also  serve  to  reduce  the  along-track  errors  as  well. 
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Figure  15:  COWPOKE  Coordinate  System  Approximation  Errors  for  HEO  Test  Case 

Despite  the  errors  present  in  the  HEO  test  case,  the  COWPOKE  equations  still  model  bulk  of  the  relative  motion 
fairly  well,  and  the  geometric  foundation  on  which  the  equations  are  built  can  provide  a  fair  amount  of  physical 
insight  into  the  relative  motion.  For  the  HEO  case,  the  physical  insight  is  very  similar  to  the  LEO  test  case  except 
we  now  observe  that  the  true  anomaly  does  not  vary  linearly  with  time.  While  the  series  expansions  do  not  provide 
a  great  deal  of  insight  into  how  the  true  anomaly  varies,  those  who  are  familiar  with  the  concept  should  be  able  to 
anticipate  the  trends.  In  Figures  10-12,  we  see  that  as  the  satellite  moves  towards  apogee  (recall  the  initial  mean 
anomaly  was  90  deg),  the  true  anomaly  rate  slows  which  pushes  the  relative  motion  extrema  to  toward  the  time  of 
perigee  passage.  It  is  interesting  to  note  that  initial  true  anomaly  difference  is  a  function  of  the  mean  anomaly 
difference,  eccentricity,  and  the  initial  mean  anomaly;  a  small  mean  anomaly  difference  at  perigee  of  a  highly 
eccentric  orbit  maps  into  a  much  larger  hue  anomaly  difference  than  if  the  initial  conditions  occurred  at  apogee. 
Additionally,  one  must  be  aware  of  this  fact  when  considering  error  sources  in  the  COWPOKE  formulation  since 
only  first  order  differences  in  the  Keplerian  elements  have  been  included  in  the  equations  above.  Finally,  from  the 
COWPOKE  equations,  one  can  see  that  the  argument  of  perigee  also  plays  a  significant  role  in  the  relative  motion 
since  the  cross-track  motion  is  dependent  on  both  the  true  anomaly  and  the  argument  of  perigee.  The  argument  of 
perigee  essentially  controls  the  phasing  between  the  cross-track  and  radial/along-track  components  of  motion. 


CONCLUSIONS 


This  paper  derived  the  Cluster  Orbits  With  Perturbations  Of  Keplerian  Elements  (COWPOKE)  equations  for 
unperturbed  satellite  motion.  A  general  framework  is  provided  to  generate  a  set  of  equations  which  describe  the 
relative  motion  between  satellites  in  eccentric  orbits  explicitly  as  a  function  of  the  Keplerian  elements  of  the 
reference  satellite,  Keplerian  element  differences,  and  time.  This  is  accomplished  using  a  geometric  description  of 
the  separation  between  the  satellites  to  form  the  basis  of  motion  and  then  by  using  a  Eourier-Bessel  series  expansion 
of  the  true  anomaly  in  terms  of  the  mean  anomaly.  True  anomaly  and  radial  differences  are  derived  using  first  order 
perturbation  methods.  The  resulting  equations  are  meant  to  provide  accurate  representations  of  the  relative  motion 
between  satellites,  to  be  simple  to  implement,  and  to  provide  physical  insight  into  the  relative  motion.  Test  cases 
show  that  this  has  been  accomplished  to  a  certain  degree. 

The  geometric  basis  on  which  the  COWPOKE  equations  are  built  is  accurate  to  the  1%  level  when  compared  to  the 
separation  distances.  This  is  currently  the  theoretical  accuracy  limitation.  For  most  applications  where  the 
eccentricities  are  below  0. 1 ,  the  COWPOKE  equations  are  easily  realizable  and  can  be  as  accurate  as  the  limitations 
of  the  geometric  basis;  however,  the  practical  limitation  for  highly  eccentric  cases  comes  from  generating  enough 
terms  in  the  series  expansion  of  the  true  anomaly  to  reduce  approximation  errors  to  an  acceptable  level.  This  work 
shows  that  the  radial  component  of  the  relative  motion  in  our  spherical  reference  system  is  particularly  sensitive  to 
true  anomaly  errors  for  high  eccentricity  orbits. 

RECOMMENDATIONS 

The  next  step  in  the  development  of  the  COWPOKE  equations  will  be  to  include  dynamic  perturbation  effects  into 
the  relative  motion.  Reference  8  has  shown  how  secular  effects  due  to  Earth’s  oblateness  can  be  incorporated  using 
Keplerian  elements  and  element  differences.  Simulations  will  also  be  performed  to  determine  what  other 
perturbations  are  required  to  support  applications  such  as  modeling  geosynchronous  clusters.  We  will  also  look  into 
the  invertibility  of  the  COWPOKE  equations.  Here  we  will  attempt  to  solve  for  Keplerian  element  differences 
based  on  desired  or  observed  relative  motion.  This  would  be  tremendously  useful  for  formation  design  and 
geosynchronous  close  approach  calculations. 
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